package de.in.tum.msspp.fail;
import java.util.ArrayList;
import java.util.List;
import java.util.SortedSet;
import java.util.TreeSet;

import Jama.EigenvalueDecomposition;
import Jama.Matrix;


public class PCAShrink {

	Matrix covMatrix;
	EigenvalueDecomposition eigenstuff;
	double[] eigenvalues;
	Matrix eigenvectors;
	SortedSet<PrincipleComponent> principleComponents;
	double[] means;
	double[][] input;
	int numComponents;

	public PCAShrink(double[][] input) {
		this.input = input;
		means = new double[input[0].length];
		double[][] cov = getCovariance(input, means);
		covMatrix = new Matrix(cov);
		eigenstuff = covMatrix.eig();
		eigenvalues = eigenstuff.getRealEigenvalues();
		eigenvectors = eigenstuff.getV();
		double[][] vectors = eigenvectors.getArray();
		numComponents = eigenvectors.getRowDimension();
		
		principleComponents = new TreeSet<PrincipleComponent>();
		for (int i = 0; i < numComponents; i++) {
			double[] eigenvector = new double[numComponents];
			for (int j = 0; j < numComponents; j++) {
				eigenvector[j] = vectors[i][j];
			}
			principleComponents.add(new PrincipleComponent(eigenvalues[i],
					eigenvector));
		}
		}

	public double[] getMeans() {
		return means;
	}
	
	public double[][] calculate(){
		Matrix originalData = new Matrix(input);
		bug("Here is the original data before adjusting by each dimension mean:");
		originalData.print(8, 4);

		//bug("There are " + numComponents + " components");
		int k = 10;
		List<PrincipleComponent> mainComponents = getDominantComponents(k);
		int counter = 1;
		//bug("Showing top " + k + " principle components.");
		for (PrincipleComponent pc : mainComponents) {
			//bug("Component " + (counter++) + ": " + pc);
		}
		Matrix features = getDominantComponentsMatrix(mainComponents);
		//bug("Feature matrix (k=" + k + ") :");
		//features.print(8, 4);

		Matrix featuresTransposed = features.transpose();
		bug("Transposed feature matrix (k=" + k + ") :");
		featuresTransposed.print(8, 4);

		double[][] matrixAdjusted = getMeanAdjusted(input, getMeans());
		Matrix adjustedInput = new Matrix(matrixAdjusted);
		//bug("Original input adjusted by dimension means (k=" + k + ") :");
		//adjustedInput.print(8, 4);
		Matrix transFormedData = featuresTransposed.times(adjustedInput.transpose());
		bug("Transformed data into PCA-space (k=" + k + ") :");
		transFormedData.transpose().print(8, 4);
		bug("*************************************************");
		
		return transFormedData.transpose().getArrayCopy();
		
	}

	/**
	 * Subtracts the mean value from each row. The means must be precomputed,
	 * which you get for free when you make a PCA instance (just call
	 * getMeans()).
	 * 
	 * @param input
	 *            Some data, where each row is a sample point, and each column
	 *            is a dimension.
	 * @param mean
	 *            The means of each dimension. This could be computed from
	 *            'input' directly, but for efficiency's sake, it should only be
	 *            done once and the result saved.
	 * @return Returns a translated matrix where each cell has been translated
	 *         by the mean value of its dimension.
	 */
	public static double[][] getMeanAdjusted(double[][] input, double[] mean) {
		int nRows = input.length;
		int nCols = input[0].length;
		double[][] ret = new double[nRows][nCols];
		for (int row = 0; row < nRows; row++) {
			for (int col = 0; col < nCols; col++) {
				ret[row][col] = input[row][col] - mean[col];
			}
		}
		return ret;
	}

	/**
	 * Returns the top n principle components in descending order of relevance.
	 */
	public List<PrincipleComponent> getDominantComponents(int n) {
		List<PrincipleComponent> retVal = new ArrayList<PrincipleComponent>();
		int count = 0;
		for (PrincipleComponent pc : principleComponents) {
			retVal.add(pc);
			count++;
			if (count >= n) {
				break;
			}
		}
		return retVal;
	}

	public static Matrix getDominantComponentsMatrix(
			List<PrincipleComponent> dom) {
		int nRows = dom.get(0).eigenVector.length;
		int nCols = dom.size();
		Matrix matrix = new Matrix(nRows, nCols);
		for (int col = 0; col < nCols; col++) {
			for (int row = 0; row < nRows; row++) {
				matrix.set(row, col, dom.get(col).eigenVector[row]);
			}
		}
		return matrix;
	}

	public int getNumComponents() {
		return eigenvalues.length;
	}

	public static class PrincipleComponent implements
			Comparable<PrincipleComponent> {
		public double eigenValue;
		public double[] eigenVector;

		public PrincipleComponent(double eigenValue, double[] eigenVector) {
			this.eigenValue = eigenValue;
			this.eigenVector = eigenVector;
		}

		public int compareTo(PrincipleComponent o) {
			int ret = 0;
			if (eigenValue > o.eigenValue) {
				ret = -1;
			} else if (eigenValue < o.eigenValue) {
				ret = 1;
			}
			return ret;
		}

		public String toString() {

			return "Principle Component, eigenvalue: " + eigenValue
					+ ", eigenvector: [" + printVector(eigenVector) + "]";

		}

		private String printVector(double[] eigenVector) {
			String vectorString = "";
			for (int i = 0; i < eigenVector.length; i++) {
				if (i != 0)
					vectorString = vectorString + ";";
				vectorString = vectorString + Double.toString(eigenVector[i]);
			}
			return vectorString;
		}
	}

	public static double[][] getCovariance(double[][] input, double[] meanValues) {
		//Zeilen
		int numDataVectors = input.length;
		//Spalten
		int n = input[0].length;

		//calculate mean and compute new matrix by it
		double[] sum = new double[n];
		double[] mean = new double[n];
		for (int i = 0; i < numDataVectors; i++) {
			double[] vec = input[i];
			for (int j = 0; j < n; j++) {
				sum[j] = sum[j] + vec[j];
			}
		}
		for (int i = 0; i < sum.length; i++) {
			mean[i] = sum[i] / numDataVectors;
		}

		//calculate C:
		//meanfree input *  meanfree input^t
		double[][] retVal = new double[n][n];
		for (int i = 0; i < n; i++) {
			for (int j = i; j < n; j++) {
				double v = getCovariance(input, i, j, mean);
				retVal[i][j] = v;
				retVal[j][i] = v;
			}
		}
		if (meanValues != null) {
			System.arraycopy(mean, 0, meanValues, 0, mean.length);
		}
		return retVal;
	}

	/**
	 * Gives covariance between vectors in an n-dimensional space. The two input
	 * arrays store values with the mean already subtracted.
	 */
	private static double getCovariance(double[][] matrix, int colA, int colB,
			double[] mean) {
		double sum = 0;
		for (int i = 0; i < matrix.length; i++) {
			double v1 = matrix[i][colA] - mean[colA];
			double v2 = matrix[i][colB] - mean[colB];
			sum = sum + (v1 * v2);
		}
		int n = matrix.length;
		double retVal = (sum / (n - 1));
		return retVal;
	}
	
	private static void bug(String what) {
		System.out.println(what);
	}
}
